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ABSTRACT 

We use full three-body scattering experiments to study the ejection of hypervelocity stars (HVSs) by 
massive black hole binaries (MBHBs) at the center of galaxies. Ambient stars drawn from a Maxwellian 
distribution unbound to the binary are expelled by the gravitational slingshot. Accurate measurements of 
thermally averaged hardening, mass ejection, and eccentricity growth rates {H, J, and K) for MBHBs in 
a fixed stellar background are obtained by numerical orbit integration from initial conditions determined 
by Monte Carlo techniques. Three-body interactions create a subpopulation of HVSs on nearly radial 
orbits, with a spatial distribution that is initially highly flattened in the inspiral plane of the MBHB, 
but becomes more isotropic with decreasing binary separation. The degree of anisotropy is smaller for 
unequal mass binaries and larger for stars with higher kick velocities. Eccentric MBHBs produce a 
more prominent tail of high-velocity stars and break planar symmetry, ejecting HVSs along a broad 
jet perpendicular to the semimajor axis. The jet two-sidedness decreases with increasing binary mass 
ratio, while the jet opening- angle increases with decreasing kick velocity and orbital separation. The 
detection of a numerous population of HVSs in the halo of the Milky Way by the next generation of large 
astrometric surveys like GAIA may provide a unique signature of the history, nature, and environment 
of the MBH at the Galactic center. 
Subject headings: black hole physics - methods: numerical - stellar dynamics 



1. INTRODUCTION 

Massive black holes (MBHs) are a ubiquitous compo- 
nent of nearby galaxy nuclei (e.g. Magorrian et al. 1998), 
and galaxies experience multiple hierarchical mergers dur- 
ing their lifetime. Following the merger of two halo-|-MBH 
systems of comparable mass ("major mergers"), dynam- 
ical friction is known to effectively drag in the satellite 
halo (and its MBH) toward the center of the more mas- 
sive progenitor: this will lead to the formation of a bound 
MBH binary (MBHB) in the violently relaxed core of 
the newly merged stellar system (Begelman, Blandford, 
& Rees 1980). Even in the case of unequal-mass merg- 
ers, gas cooling appears to facilitate the pairing process 
by increasing the resilience of the companion galaxy to 
tidal disruption (Kazantzidis et al. 2005). It is expected 
then that many galaxies will host wide MBHBs during 
cosmic history (e.g. Volonteri, Haardt, & Madau 2003). 
As the binary separation decays, the effectiveness of dy- 
namical friction slowly declines because distant stars per- 
turb the binary's center of mass but not its semi-major 
axis. The bound pair then loses orbital energy by cap- 
turing stars passing in its immediate vicinity and ejecting 
them at much higher velocities (gravitational slingshot). 

It was first pointed out by Hills (1988) that the tidal 
breakup of binary stars by a MBH at the Galactic cen- 
ter may eject one member of the binary with velocities 
^ 1000 km s~^. Such 'hypervelocity stars" (HVSs) are 
also produced by three-body interactions between ambi- 
ent stars with low angular momentum orbits and a "hard" 
MBHB, i.e. a binary whose binding energy per unit mass 
exceeds the star specific kinetic energy. Assuming SgrA* 
to be one component of a MBHB, Yu & Tremaine (2003) 



estimated the number of HVSs expected within the solar 
radius to be ^ 10^. Brown et al. (2005) reported the first 
discovery of a HVS (with a Galactic rest-frame velocity in 
excess of 700 kms~^) in the Galactic halo. This and five 
more HVSs, recently discovered by Hirsch et al. (2005) 
and Brown et al. (2006a, 2006b), are all consistent with a 
Galactic center origin, while an ejection from the LMC is 
more plausible for the seventh HVS currently known (Edel- 
mann et al. 2005). HoUey-Bockelmann et al. (2005) have 
proposed that the anomalously fast intracluster planetary 
nebulae identified in the Virgo Cluster by Arnaboldi et al. 
(2004) may also be associated with close three-body inter- 
actions with a MBHB at the center of M87. Unbound 
HVSs travel with velocities so extreme that dynamical 
ejection from a relativistic potential is the most plausi- 
ble origin, and are becoming increasingly recognized as an 
important tool for understanding the history, nature, and 
environment of nuclear MBHs. 

Stars expelled by a MBHB are expected to form a sub- 
population with very distinct kinematics (e.g. Quinlan 
1996, hereafter Q96) as well as spatial structure (Zier & 
Biermann 2001). The phase-space distribution of HVSs 
ejected by an intermediate-mass black hole (IMBH) inspi- 
ralling into SgrA* has been recently studied analytically 
by Levin (2005). Stars bound to SgrA* and drawn from an 
isotropic cusp are ejected in a burst lasting a few dynami- 
cal friction timescales: most stars are expelled isotropically 
if the inspiral is circular, or in a broad "jet" aligned with 
the IMBH velocity at pericenter if the inspiral is eccen- 
tric (Levin 2005). Yet most HVSs are produced during 
the phases of the inspiral that cannot be modeled analyt- 
ically. This is the first paper in a series aimed at a de- 
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tailed numerical study of the interaction of MBHBs with 
their dense stellar environment. We use here full three- 
body scattering experiments of the ejection of background 
stars by MBHBs at the center of galaxies to address the 
kinematic properties of HVSs. Ambient stars drawn from 
a Maxwellian distribution unbound to the binary are ex- 
pelled by the gravitational slingshot. Reaction rates are 
obtained by numerical orbit integration from initial con- 
ditions determined by Monte Carlo techniques. The plan 
of the paper is as follows. In § 2 we describe our suite of 
three-body scattering experiments. In § 3 we present accu- 
rate measurements of the binary hardening, mass ejection, 
and eccentricity growth rates for MBHBs embedded in a 
fixed stellar background, reproducing Q96's classical re- 
sults. In § 4 we discuss the detailed kinematic properties 
of the ejected subpopulation. Finally, we present our con- 
clusions in § 5. 

2. SCATTERING EXPERIMENTS 

Consider a binary of mass M = Mi + M2 (M2 < Mi), 
reduced mass fi — M1M2/M, and semimajor axis a, orbit- 
ing in the (x, y) plane in a background of stars of mass m*. 
In the case of a light intruder with m* <C M2, the problem 
is greatly simplified by setting the center of mass of the 
binary at rest at the origin of the coordinate system. It is 
then convenient to define an approximated dimensionless 
energy change C and angular momentum change B in a, 
single binary-star interaction as (Hills 1983) 



C = 



and 



B 



M 


AE 


2to* 


E 


M AL, 



L, 



aAE^ 



M AL^ 



(1) 



(2) 



Here AE/E is the fractional increase (decrease if negative) 
in the orbital specific binding energy E = —GM/{2a), 
ALz/Lz is the fractional change in orbital specific an- 
gular momentum Lz = \jGMa{l — e^), while A_B* and 
ALz* are the corresponding changes for the interacting 
star. Conservation of total energy and angular momentum 
lead to the following expression for the change in orbital 
eccentricity e. 
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The quantities B and C are of order unity and are de- 
rived by three-body scattering experiments that treat the 
star-binary encounters one at a time (Hut & Bahcall 1983; 
Q96). For each encounter one solves nine coupled, second- 
order, differential equations. 
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supplied by 18 initial conditions Yi{t = 0), r'j(t = 0). The 
incoming star is moved from r — 00 to r^ = [10^°/i/M]^'*a 
on a Keplerian orbit about a point mass M. At r^ the force 
induced by the quadrupole moment of the binary is 10 or- 
ders of magnitude smaller than the total force acting on 

^See http;//www.unige.ch/math/folks/hairer /software. html. 



the star at a distance a, and numerical integration starts. 
The initial conditions define a point in a nine-dimensional 
parameter space given by: 

• the mass ratio q ~ M2/M1 of the binary; 

• the eccentricity of the binary e; 

• the mass of the star ?n, ; 

• the asymptotic initial speed of the incoming field star 



• the impact parameter at infinity b (the distance at 
which the star would pass the binary if it fell no 
attraction); 

• four angles: 9 and (j) describing the initial direction 
of the impact, ■0 its initial orientation, and ^ the 
initial binary phase. 

A significant star-binary energy exchange (i.e. character- 
ized by a dimensionless energy change C > 1) occurs only 
for v/Vc < A/M2/M, where V^ = ^GMja is the binary 
orbital velocity (the relative velocity of the two holes if the 
binary is circular, see, e.g., Saslaw, Valtonen, & Aarseth 
1974; Mikkola & Valtonen 1992). We sample this quantity 
in the range 3 x lO^'V^a/M < w/K < 30^M2/M with 
80 logarithmically equally spaced grid points. Note that 
varying the incoming velocity at a fixed Vc is equivalent 
to varying the binary separation a at a fixed v. When av- 
eraged over a Maxwellian distribution, the range of v jVc 
sampled in our experiments therefore probes the slingshot 
mechanism over about four decades in binary separation. 
Each scattering experiment requires five random num- 
bers. The impact parameter is related to the star pericen- 
ter distance r^ by 



2GM\ 



(5) 



For each incoming velocity the impact parameter is ran- 
domly sampled according to an equal probability distribu- 
tion in 6^ (equivalent to a probability weight proportional 
to 6), in an interval corresponding to a range in scaled 
pericenter distance r^ja of [0, 5]. In a number of test cases 
we have extended this range to [0, 10], finding little differ- 
ences in the measured hardening rates. The velocity angles 
Q and are randomly generated to reproduce a uniform 
density distribution over a spherical surface centered in the 
origin of coordinate system, while the orientation angle '0 
is chosen from a uniform distribution in the range [0, 27r]. 
The probability distribution of the binary phase ^ is sam- 
pled by weighting randomly selected angles according to 
the time spent by the binary in any given phase range. 
Given a set of initial conditions, we integrate the system 
of differential equations (4) using the most recent version 
of the subroutine DOPRK'^, which is based on an explicit 
Runge-Kutta method of order 4(5) due to Dormand & 
Prince (1978). A complete description of the integrator 
can be found in Hairer, Norsett, & Wanner (1993). 

As the total energy of the binary-star system is always 
negative, a bound triple system may form temporarily, but 



will typically dissolve within 10 to 100 crossing times (Val- 
tonen & Aarseth 1977). The integration is stopped if one of 
the following events occurs: (a) the star leaves the sphere 
of radius r^ with positive total energy; (b) the integra- 
tion reaches 10® time steps (corresponding to between a 
few hundred and a thousand binary orbital periods); or 
(c) the physical integration time exceeds 10"'^'' yrs. The 
integration of the full three-body problem allows us to di- 
rectly control the conservation of total energy and angular 
momentum. The code adjusts the integration stepsize to 
keep the fractional error per step, e, in the position and 
velocity below a level which was set to 10^^^. This allows 
a total energy conservation accuracy AE/ E ^ 10~^ in a 
single orbit integration, while, for m^,/M ~ 10"^, the star 
energy is conserved, in a single orbit, at level of one part 
in a hundred. We have also checked that our choice of 
erasing from the record stars that get captured for long 
times does not affect results appreciably. Abandoned inte- 
grations involve, for the major part, encounters with stars 
that get captured in very weakly bound orbits (energet- 
ically "poor" events) and make many revolutions before 
being expelled. The fraction of erased stars that are in- 
stead captured in tightly bound orbits (hence energetically 
"rich") is quite small, <^ 10~^, assuring that the global im- 
pact of such encounters on the evolution of binary orbital 
parameters is indeed negligible. 

We have performed 24 sets of scattering experiments 
for binary mass ratios q^l, 1/3, 1/9, 1/27, 1/81, 1/243, 
and initial eccentricities e=0.01, 0.3, 0.6, 0.9. After each 
orbit integration the x, y, z components of v and L are 
stored. Each run involves 4 x 10® stars: we collect par- 
tial outputs after 5 x 10^ orbit integrations correspond- 
ing to a given initial speed u, calculate the values of (C) 
and {B) averaged over orbital angular variables, and then 
evaluate the hardening, eccentricity growth, and mass ejec- 
tion rates as defined in the § 3. Statistical errors are es- 
timated by evaluating the rates from ten different orbit 
subsets, and then computing standard deviations (see § 3). 



Fig. 1. — Top panel: mean energy exchange (C) as a function of 
dimensionless impact parameter x for an equal mass binary with ec- 
centricity e = (solid line), 0.3 (long-dashed line), 0.6 [short-dashed 
line), and 0.9 [dot-dashed line). Bottom panel: same, but for a cir- 
cular binary with mass ratio q = 1/3 (solid line), 1/9 [long-dashed 
line), 1/27 (short-dashed line), and 1/81 (dot-dashed line). 
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Fig. 2. — Mean angular momentum minus energy exchange 
(B — C) parameter as a function of dimensionless impact parameter 
X for a binary with eccentricity e = (solid line), 0.3 [long-dashed 
line), 0.6 (short-dashed line), and 0.9 [dot-dashed line). Top panel: 
equal-mass binary [q = 1). Bottom panel: q = 1/9. 

Figure 1 shows examples of the energy exchange param- 
eter (C) versus the star impact parameter h in the limit 
V <^Vc i.e. for a "hard" binary. We have defined a rescaled 
dimensionless impact parameter x as in Q96: 



= b/ha, bl = 2GMa/v^ 



(6) 




where 6o is the impact parameter at infinity leading to 
a star pericenter rp = a including gravitational focusing. 
For a circular, equal-mass binary, (C) has a maximum 
for impact parameter equal to the distance of each binary 
member from the center of mass (x ~ 0.5), and then de- 
creases reaching a constant value (C) ~ 1 for small im- 
pact parameters. In the case of highly eccentric binaries 
(e — 0.9), the energy exchange is almost monotonic: a star 
with small impact parameter interacts with the pair near 
pericenter when it has maximum speed, and this results in 
a large energy exchange. For circular binaries and x <^ 1, 
one finds 1.5 <; (C) <; 2 nearly independently of the mass 
ratio. The peak at a; = 1 for g <C 1 corresponds to the 
location of M2, as Mi is now at rest and M2 orbits around 
Ml at distance a. From equation (1) it follows that stars 
typically experience an energy gain proportional to the re- 
duced mass of the binary, and inversely proportional to 
the binary separation. The combination {B{x) — C{x)), 
which governs the eccentricity evolution (see eq. 12), is 
shown in Figure 2. For a circular binary {B — C) is zero 
independently of impact parameter, i.e. circular binaries 
remain circular. If the binary is eccentric, (B — C) shows a 
positive peak for impact parameters corresponding to close 
encounters with M2 or with one member of an equal-mass 



Massive black hole binaries 



pair. Interactions with 0.5 ^ a; ^ 2 stars tend to increase 
the eccentricity of the binary, whereas stars with low im- 
pact parameter have either a small effect (q — 1/9) or can 
significantly decrease the eccentricity (q — 1). Typically, 
B ex ALz* assumes positive values, i.e. scattered stars 
gain angular momentum along the z-axis and will tend to 
corotate with the binary. 

3. HARDENING IN A FIXED BACKGROUND 

As described in Q96, the binary evolution in an isotropic 
fixed background of stars of density p and one-dimensional 
velocity dispersion a at infinity is determined by three di- 
mensionless quantities: the hardening rate 



Gpdt \aj ' 



(7) 



the mass ejection rate (Mcj is the stellar mass ejected by 
the binary) 

Mdln(l/a)' ^ ' 

and the eccentricity growth rate 

de 



K = 
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The average hardening rate for a Maxwellian stellar veloc- 
ity distribution f{v,a) = (27rcr^)^'^/^ exp(— i;^/2cr^) is 



where 



H{(j)= f{v,(7)-Hl{v)4:TTV^dv, (10) 



(11) 



Hi{v) = 87r / {C)xdx 



is the dimensionless hardening rate if all stars have the 
same velocity v. An expression analogous to equation 
(10) relates the thermally-averaged eccentricity growth 
rate K{<t) to Ki{v), where 
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Both Hi and Ki are independent of M and tti* . Figure 3 
shows the hardening rate H versus binary separation a 
as derived from our scattering experiments. As found by 
Q96, H is approximatively constant for 



a < Gh 



GM2 
1^ 



(13) 



impact parameters at increasing eccentricity. Note that 
H drops to zero for a ^ lOa/j. We choose the follow- 
ing ejection criterion to measure the rate at which the 
binary ejects stars. The binary is assumed to be em- 
bedded in a bulge of mass Mb and stellar density pro- 
file approximated by a singular isothermal sphere (SIS). 
Stars arc counted as "ejected" from the bulge if, after 
three-body scattering, their velocity V far away from the 
binary is greater than the escape velocity from the ra- 
dius of influence of the binary, Hnf = GM/ (2(T^). The 



SIS potential is 



^2ct2 [\n{GMB/2a^ 



1] (for 



r < Re 
is then 



GMs/'ia ), and the escape speed from n, 



v/^ = 2ay/[\n{MB/M) + 1] = 5.5a, (14) 



where the second equality comes from the adopted bulge- 
black hole mass relation M = 0.0014 Mb (Haring & Rix 
2004). Note that this yields a more conservative ejection 
criterion than the conventional choice V > Uosc = 2\/6a 
adopted by Q96. ^ Denoting with Fcj{x,v,a) the frac- 
tion of incident stars with impact parameter x and initial 
velocity v that satisfy equation (14) after a three-body in- 
teraction, the thermally averaged ejection rate can then 
be written as 



J(a) 



1 
H 



f{v, cr) — Attv'^ dv 47r 



Fcjxdx. (15) 



We have found that the (Maxwellian averaged) rates H, J, 
and K derived from scattering experiments can be fitted 
to within few percent by the following functions: 



J 



and 



H = A{l+a/ao)\ 
K = A{l+a/aoy +B. 



(16) 
(17) 

(18) 



The parameters of the fits to H and J are listed in Tables 
1 and 2 for different binary mass ratios and circular orbits, 
while fit parameters to K are listed in Table 3 for differ- 
ent values of q and e. The binary eccentricity growth and 
mass ejection rates K and J (averaged over a Maxwellian 
distribution) are plotted in Figures 4 and 5 as a function 
of a/afi- The parameter K is close to zero for a ^ ah, and 
grows monotonically in the case of eccentric orbits as the 
binary shrinks. Smaller values of q lead to larger growth 
rates. In the circular case K is negligible at every sep- 
aration, i.e. circular binaries remain circular. The mass 
ejection rate has only a weak dependence on eccentricity 
and mass ratio. 



Defining as "hard" a binary whose orbital separation is 
smaller than ah, it follows then that a "hard binary hard- 
ens at a constant rate". Note that there is no explicit 
dependence on a once the rate is expressed as a function 
of a I ah- H is found to be a decreasing function of q but 
increases with increasing eccentricity. The latter trend 
can be understood from Figure 1, which shows a signi- 
ficative enhancement in the energy change (C) at small 

^When the binary first becomes hard, only a few stars acquire a kick velocity large enough to escape the host bulge. Many scattered 
stars will instead return to the central region on nearly unperturbed, small impact parameter orbits, and will undergo a second super-elastic 
scattering with the binary. We will quantify the role of these "secondary slingshots" in determining the hardening of the pair in a subsequent 
work. 
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Table 1 

Best fit parameters describing (sec eq. 16) the hardening rate H 

for a circular binary with varying mass ratio q. The parameter ao 

is given in units of a^. 
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Table 2 

Best fit parameters describing (see eq. 17) the mass ejection rate J 

for a circular binary with varying mass ratio. The parameter ao is 

given in units of a^ . 
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Table 3 

Best fit parameters describing (see eq. 18) the eccentricity growth 

rate K for varying initial eccentricity e and mass ratio q. The 

parameter ao is given in units of a^. 

It is interesting to compare our fits to those provided 
by Q96. While our results for the dimensionless rates 
Hi and Ki are consistent, within the statistical errors, to 
those of Q96, we have fit directly the Maxwellian- averaged 
rates for ease of use. Our choice of a different velocity 
threshold for ejection makes our mass ejection rate J(cr) 
generally lower than that of Q96. Regarding H{a), lim- 
iting values for small separations are virtually identical 
to those of Q96, while at a ~ a/j the fit in Q96 is 20% 
above our formula. As for the eccentricity growth rate, we 
have covered a much wider region of the (q, e) parameter 
space, compared to the limited sampling provided by Q96. 



Fig. 3. — Binary hardening rate H (averaged over a Maxwellian 
velocity distribution) versus a/a/^. Left panel: e =0, 0.3, 0.6, 0.9 for 
5 = 1. Right panel: q =1/3, 1/9, 1/27, 1/81 for e = 0. Line style as 
in Fig. 1. 
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Fig. 4. — Binary eccentricity growth rates as a function of bi- 
nary separation for different mass ratios. Different line styles are 
for e =0, 0.3, 0.6, 0.9, as in the top panel of Fig. 1. 



4. KINEMATICS OF HYPERVELOCITY STARS 

When the MBHB separation is a <; ah, only a small frac- 
tion (<J 1%) of bulge stars have low-angular momentum 
trajectories with pericenters lying within a (the binary's 
geometrical "loss cone"). If the loss-cone is constantly re- 
filled as the pair shrinks, then a substantial subpopulation 
of suprathermal HVSs will be produced via the gravita- 
tional slingshot. In this section we use our scattering ex- 
periments to study the kinematic of HVSs as a function 
of binary mass ratio, eccentricity, and orbital separation. 
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where w = V/Vc, and h = ^/2q/{l + q). Best fitting pa- 
rameters are given in in Table 4 for different values of the 
eccentricity, while the fitting formula and the scattering 
experiment data are compared in Figure 8. The velocity 
distribution of HVSs {V > Wosc) at a given binary separa- 
tion is related to the mass ejection rate by 



J{a) 



d{M,-JM) 
da/ a 



o/l4 



f{w)dw. 



(20) 



Equation (20) sets the normalization constant A in equa- 
tion (19). 



Fig. 5. — Mass ejection rate J versus a/a/^. Top panel: e =0, 0.3, 
0.6, 0.9 for q = l. Bottom panel: q =1/3, 1/9, 1/27, 1/81 for e = 0. 
Line style as in Fig. 1. 



4.1. Velocity distribution 

In a scattering event, a star that starts with a low 
initial velocity v passes the two MBHs at a distance 
~ a and leaves with a gain to its kinetic energy, V"^ ~ 
2C{iJ,lM)V^ = [8CcrV(l + <7)](a,i/a) (eq. 1), where C de- 
pends on the impact parameter (see Fig. 1). Figure 6 
shows the average final velocity of scattered stars (V) as 
a function of inverse binary separation, for different mass 
ratios and eccentricities. The population includes all stars 
with maximum impact parameter at infinity x < 2, cor- 
responding to a pericenter rp < 4a. The curves clearly 

follow the expected V ex ^Jofja scaling. For Ohja ^ 10, 
most scattering events produce HVSs that escape from the 
bulge. Neither the binary mass ratio nor its eccentricity 
have a large effect on the average final velocity. Larger 
eccentricities produce a more prominent high- velocity tail 
of HVSs as in this case the orbital velocity of M^ close to 
pericenter is larger than Vc, allowing for more energetic 
slingshots. Our numerical experiments show that both a 
small mass ratio and a large eccentricity tend to increase 
the fraction of extremely energetic scattering events. 

The properties of scattered stars are best described 
noticing that, when stellar ejection velocities are measured 
in units of binary orbital velocity Vc^ the high- velocity 
tail of the resulting distribution function is actually scale- 
invariant, i.e. independent of binary separation for fixed 
q and e. This is clearly shown in Figure 7. Scale in- 
variance is broken by the choice of an absolute velocity 
threshold (e.g. V > Ucsc) that, in units of Vc, is a de- 
creasing function of Oh/a. The high-velocity tail of the 
distribution function of scattered stars can be described 
as a broken power-law in the range Acr <^ V <^ SVmax, 

where Knax = K ^(1 + e)/(l — e)/(l + q) is the velocity 
of the lighter black hole at the pericenter: 

f{w)^^{w/hr[l + {w/hY']\ (19) 
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f^ 
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0.236 


-0.917 


16.365 


-0.165 


0.3 


0.242 


-1.067 


11.722 


-0.235 


0.6 


0.385 


-0.765 


4.627 


-0.726 


0.9 


0.556 


-0.599 


2.375 


-1.420 



Table 4 

Best fit parameters describing (see eq. 19) the velocity distribution 
of scattered stars with 4cr < V < 3ymax. 
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Fig. 6. — Average final velocity (V) of all scattered stars versus 
inverse binary separation ai^/a. Top panel: eccentric e = 0.9 binary 
with g = 1 {solid lines), q = 1/9 {long-dashed lines), and q = 1/81 
{short-dashed lines). The lower set of curves shows results for the 
entire population of scattered stars (defined as those with impact 
parameter at infinity x < 2), while the upper set of curves shows 
only HVSs with V > i)csc = 5.5(7. The dotted line in both panels 
marks the value V = Vosc. Bottom panel: same but for q = 1/9 
binary with e = {solid lines), e = 0.3 {long-dashed lines), e = 0.6 
{short-dashed lines), and e = 0.9 {dot-dashed lines). 
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Fig. 7. — Differential distribution of V/Vc for all stars scattered at 
binary separation ay^ja = 1 (solid line), 10 (long-dashed line), and 
100 (short- dashed line). As the binary shrinks the peak shifts to 
increasingly small values of V/Vc as injected stars are drawn from 
a thermal distribution of fixed velocity dispersion a. Left panel: bi- 
nary with q = 1 and e = 0.9. Right panel: binary with q = 1/81 
and e = 0.9. 
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Fig. 8. — Differential distribution of V/Vc for all stars scattered 
at binary separation a/^/a = 100. In each panel values are shown, 
right to left, for binary mass ratios q = 1, 1/9 and 1/81. Scattering 
experiment data are shown as solid lines, while fitting formula re- 
sults, extended up to V = SVmax (eq- 19, and Table 4), as dashed 
lines. The appropriate values of the eccentricity are labeled in each 
panel. 



4.2. Angular distribution 

Analytic expressions for the time-dependent phase- 
space distribution of stars ejected from the Galactic cen- 
ter as a result of inspiral of an intermediate-mass black 
hole {q <C 1) have been recently derived by Levin (2005). 
HVSs are found to follow a characteristic angular pattern: 
i) they are ejected preferentially in the orbital plane of 
the MBHB, and ii) perpendicular to the semimajor axis 



of the pair in the case of eccentric binaries. For a given 
mass ratio and eccentricity, the magnitude of these three 
effects is a function of orbital separation and stellar final 
velocity. Energetic three-body interactions eject stars in 
the direction of maximum black hole orbital speed, gen- 
erating a larger anisotropy in higher velocity stars, and a 
positive z— component of the stellar angular momentum. 
For small binary separations the degree of anisotropy of 
the ejected stars is reduced as, in order to generate a kick 
above a given V , the three-body interaction needs not to 
be as strong as in the case of large separation. 

4.2.1. Latitude 

Figure 9 shows (6*^) as a function of binary separation, 
where — 7r/2 < 6 < 7r/2 is the latitude of ejected stars, 
i.e. the angle between the star velocity vector at infinity 
and the binary orbital plane. As an isotropic distribution 
would yield (0^) = 7r^/4 — 2 ~ 0.467, values lower than 
this indicates that the ejected stars tend to be flattened 
towards the binary orbital plane. We plot results for the 
population of HVSs as a whole {V > Wcsc) and for a subset 
with V > 9wesc- Note that a value (6*^) = at large separa- 
tions simply means that no stars are ejected with a velocity 
exceeding the given threshold, not that stars are actually 
scattered exactly in the binary orbital plane. Higher ve- 
locity stars are more flattened into the orbital plane, and 
smaller binary separations lead to a more isotropic angular 
distributions. At a fixed mass ratio and orbital separation, 
more eccentric binaries eject stars in a more isotropic fash- 
ion, while at a fixed eccentricity binaries with smaller mass 
ratios produces a more isotropic distribution of HVSs. 
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Fig. 9. — Mean value of 9^ versus a/^/a for HVSs, where 9 is the 
latitude angle between the star velocity vector at infinity and the 
binary orbital plane. In each panel, the upper set of curves refers 
to all ejected stars (V > t'esc), the lower set of curves to stars with 
V > 9fesc- Top panel: binary with mass ratio q = 1/9 and ec- 
centricity e = (solid lines), e = 0.3 (long-dashed lines), e = 0.6 
(short-dashed lines), and e = 0.9 (dot-dashed lines). Bottom panel: 
eccentric binary with e = 0.9 mass ratio 5 = 1 (solid lines), q = 1/9, 
(long-dashed lines), and q = 1/81 (short-dashed lines). 



4.2.2. Longitude 
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Scattered stars typically receive a kick along the direc- 
tion of maximum velocity of the MBHs. This implies that, 
in the case of eccentric binaries, the spatial distribution of 
HVSs will form a broad jet perpendicular to the semima- 
jor axis. For small binary mass ratios, the jet tends to be 
one-sided since Mi is practically at rest and the interaction 
takes place close to the pericenter of M2 (Levin 2005). For 
equal-mass binaries, we expect a two-sided symmetric jet 
instead. To quantify these effects we have assumed that 
the semimajor axis coincides with the x— axis, and that the 
two holes orbit counterclock-wise in the xy-plane. The M2 
(Ml) velocity vector at the pericenter then forms an angle 
37r/2 {tt/2) relative to the x— axis. Let (p be the star longi- 
tude, i.e. the angle between the a;— axis and the projection 
onto the orbital plane of the star velocity vector at infin- 
ity. In the case of an azimuthally symmetric distribution, 
{(j)) ~ TT, with relative dispersion a^/{(f>) = l/-\/3 = 0.58. 
Figure 10 shows the mean longitude of HVSs and its dis- 
persion for an eccentric binary with different mass ratios, 
as a function of a^ja^ both for the whole population of 
HVSs stars {V > Wosc) and for a a subset with V > 9fcsc. 
In the case g <C 1 the jet is oriented along the velocity 
of M2 at pericenter, while for q = 1 the value {(j)) ~ w 
does not denote axisymmetry but rather a two-sided jet 
oriented perpendicular to the semimajor axis (the bottom 
panel shows that the angular dispersion is lower than its 
symmetric value). HVSs are ejected more symmetrically 
as the binary shrinks and, at any given separations, the 
degree of axisymmetry increases with decreasing kick ve- 
locity. The mean longitude angle as a function of eccen- 
tricity is shown in Figure 11 for a MBHB binary with 
q = 1/9. While the population of HVSs as a whole is 
nearly symmetric for all eccentricities regardless of binary 
separation, the high velocity subsample shows a clear az- 
imuthal asymmetry for ah /a ^ 20 that decreases with de- 
creasing separations. MBHBs already significantly beam 
HVSs for e — 0.3: larger eccentricities produce similar 
values of {(p) but at larger separations. 

We note that, as well as ejection velocity, also the angu- 
lar properties of scattered stars depend on the ratio V/Vc 
but not on the binary hardness ah/a. This can be clearly 
seen in Figure 12, where the mean latitude and longitude 
angles are plotted versus V/Vc for a binary with q — 1/9 
and e = 0.6, at separations a/ah — 1, 0.1 and 0.01. 
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Fig. 10. — Top panel: mean longitude (0) of ejected stars as a 
function of orbital separation. The binary semimajor axis lies along 
the x— vector, the MBHs orbit counterclock-wise in the xy— plane, 
and the binary eccentricity is set to 0.9. The solid, short-dashed, 
and long-dashed curves show values for q = 1, 1/9, and 1/81, re- 
spectively. Thin lines: entire population of ejected stars {V > ticsc)- 
Thick lines: subset of HVSs with V > 9i'osc- Bottom panel: relative 
azimuthal dispersion a^j {ip). 
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Fig. 11. — Same as Fig. 10 but for a binary with q = 1/9 and 
eccentricity e = {solid lines), 0.3 {long-dashed lines), 0.6 {short- 
dashed lines), and 0.9 {dot-dashed lines). Thin lines: entire popu- 
lation of ejected stars {V > Dose)- Thick lines: subset of HVSs with 

V > 9l)esc. 
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circular binary by Yu & Tremaine (2003) using a vari- 
ety of observational and theoretical arguments. Figure 
13 shows the quantity ci(7VHVs)/<^lii(a), the number of 
HVSs (assuming jti* = 1 Mq) ejected per logarithmic bi- 
nary separation with V > 1.68 Ucsc- The adopted velocity 
threshold implies V > 840 kms~^ at the influence radius 
of the binary: it corresponds to the escape velocity from 
the center of the MW galaxy,^ and translates within the 
MW potential to F > 450 km/s 10 kpc away from SgrA*. 



Fig. 12. — Mean values of 9^ (right scale axis) and <f> (left scale 
axis) versus V/Vc at different binary separations: a/a/^ = 1 (solid 
lines), a/a/^ = 0.1 (long-dashed lines), and a/a/^ = 0.01 (short- 
dashed lines). A mass ratio q = 1/9 and an eccentricity e = 0.6 are 
assumed. 



5. CONCLUSIONS 

We have performed full three-body scattering experi- 
ments in order to study the detailed kinematic properties 
of hypervelocity stars by massive black hole binaries at 
the center of galaxies. Ambient stars are drawn from a 
Maxwellian distribution unbound to the binary, and are 
expelled by the gravitational slingshot. Numerical orbit 
integration from initial conditions determined by Monte 
Carlo techniques provides accurate measurements of ther- 
mally averaged hardening, mass ejection, and eccentricity 
growth rates for MBHBs in a flxed stellar background. 
We have shown that binary-star interactions create a sub- 
population of HVSs on nearly radial orbits, with a spatial 
distribution that is initially highly flattened in the inspi- 
ral plane of the MBHB, but becomes more isotropic with 
decreasing binary separation. The degree of anisotropy 
is smaller for unequal mass binaries and larger for stars 
with higher kick velocities. Eccentric MBHBs produce a 
more prominent tail of high-velocity stars and break ax- 
isymmetry, ejecting HVSs along a broad jet perpendicular 
to the semimajor axis. The jet two-sidedness decreases 
with increasing binary mass ratio, while the jet opening- 
angle increases with decreasing kick velocity and orbital 
separation. 

It is interesting to quantify the properties of the HVSs 
that would populate the halo of the Milky- Way (MW) 
in the presence of a MBHB at the Galactic center. We 
assume that the binary total mass is Mi -\- M2 == 3.5 x 
10^ Mq, that Sgr A* is the most massive component Mi 
of the pair, that the binary mass ratio is q = 1/81, and 
that the loss-cone is always full. From our discussion in 
§4.1, we expect the small mass ratio to result in a large 
number of HVSs. The pair is allowed to shrink from 
a = ah down to a = 0.1 a;i due to three-body interac- 
tions, within the allowed parameter space derived for a 

^We modeled the luminous component of the MW as in Miyamoto 
(2005). 




Fig. 13. — Number of HVSs escaping the MW. The binary 
total mass is M\ + M2 = 3.5 X 10^ Mq, and the mass ratio is 
q = M2/M1 = 1/81. Solid line: binary with initial eccentricity 
e = 0. Dot-dashed line: binary with initial eccentricity e = 0.9. 
The total number of HVSs ejected as the pair shrinks from a = af^ 
to a = ah/2 (a = a^/lO) is 5 X 10* (lO"). 

The curves show results for different initial eccentrici- 
ties: in a self-consistent treatment of the evolution of a 
MBHB, we have used the results of our scattering experi- 
ments to account for the changing binary eccentricity as its 
orbit decays. As the orbital separation decreases from ah 
to aft,/2, we flnd that the total number of HVSs is 5 x 10**, 
for a total ejected mass of some 1.2 M2, independent of 
eccentricity: in the case of a circular binary the angular 
distribution of such HVSs is characterized by {9"^) — 0.3 
and a,),/ {(p) = 0.58, while for an eccentric pair with e — 0.9 
the latitude and azimuthal dispersions are {9"^) — 0.37 and 
(70 /((/)) = 0.6. In the case in which the binary separation 
shrinks from the hardening radius down to a = 0.1 o/j in- 
stead, the total number of ejected HVSs increases by a 
factor of 20, to 10®: a circular binary produces an angular 
distribution with {9"^) = 0.36 and a^/{(t)) = 0.58, while 
for an eccentric pair with e = 0.9 we flnd (0^) = 0.4 and 
a^l {(j)) ~ 0.6. The detection of a numerous population of 
HVSs in the halo of the Milky Way by the next generation 
of large astrometric surveys like GAIA may thus provide 
a clear signature of the history, nature, and environment 
of the MBH at the Galactic center. 
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